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CHAPTER I 



INTRODUCTION . 

Many types of structures are designed for use below and/or 
above the surface or on the shoreline of oceans and various large 
inland water areas. They can be floating and/or structurally supported 
to withstand gravitational forces. These structures are constructed to 
withstand wave forces to the extent present day knowledge of water 
wave mechanics, standard of practice, and appropriate factor of safety 
allow. A criterion used in design may be the "design wave", a 
statistical culmination of known wave heights recorded in the area 
combined with an appropriate crest length which is used to compute the 
largest horizontal forces which the structure will be required to 
withstand. 

An important design assumption used for many structures built 
to withstand water wave forces is the reflection coefficient. Defined 
as the ratio of the height of the wave reflected from the structure to 
the height of the incident wave coming in contact with the structure 
{1} , it is an indicator of wave energy reflected from the structure 
as the square root of wave energy is directly proportional to wave 
height or amplitude. This reflected energy is an indicator of the 
direct effect the structure will have on any neighboring structure or 
shoreline. Combined with knowledge of the transmitted energy which 
passes through the structure and possible motion of the structure 
itself, it can be used to determine energy absorbed by the structure 
and forces which the structure will be required to withstand to 
perform properly. 

A 

Numbers shown in brackets refer to the References listed on pp. 39-40. 
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Specifically, breakwaters, floating bridges, piers, wharves, 
and other structures designed to withstand and even modify water 
wave forces have a particular requirement for knowledge of reflection 
coefficients both in their design and for follow-up monitoring of true 
effects after construction is completed. Determination of reflection 
coefficients either from model study during design or from study of 
the finished structure itself is not always simple, even if the 
situation allows measurement of both the incident wave sea state 
alone and the combined reflected and incident wave sea state. These 
measurements are highly dependent on angle of attack of the incident 
wave train, local conditions (water depth, etc.), and the particular 
frequencies and amplitudes associated with different incident wave 
trains found at different times. The problem is particularly 
difficult when one is working with a structure such as a floating 
bridge where it is impossible to measure the incident wave sea state 
prior to its combination with reflected wave components. 

Morden {2, 3}, Thornton and Calhoun {4} and others have 
developed deterministic approaches to computing reflection co- 
efficients from the combined incident/reflected sea state based on 
linear deep water wave theory. Their methods require the use of two 
or more wave measurement devices in the combined wave field. 

In contrast, the thrust of the current investigation is 
the separation of incident and reflected wave energy from the combined 
sea state using only one wave measurement device and cepstral analysis 
techniques. The present use of cepstral analysis is primarily in 
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seismic epicenter estimation, speech parameter estimation, photo- 
graphic image enhancement, radar, sonar, and other fields where one 
desires separation or identification of the reflected signal from the 
original signal or the time of delay between the two signals. The 
literature contains several references to recent cepstral analysis 
work performed in the speech and seismic arenas but apparently no one 
has ever attempted application of this process to water waves, 

Cepstral analysis wave recovery and decomposition techniques 
assume neither the reflection coefficient nor the delay time are 
functions of frequency. As this is not an absolutely correct 
assumption with water waves, it is probably understandable why this 
technique has not been actually attempted with water waves. However, 
cepstral analysis decomposition would allow determination of 
reflection coefficients using only one wave measurement device 
collecting data from the combined incident/reflected wave field. It, 
therefore, would be of great benefit if this method could be reliably 
applied with respect to ocean waves. 

Chapter Two discusses very briefly some of the general 
mathematical concepts employed in this work and directs the reader to 
other references if one whishes to delve into discrete time series 
analysis in more detail. It by no means is an exhaustive explanation 
of background material but is instead a hint at what the field 
involves . 

Chapter Three contains a description of cepstral analysis 
in general. Starting with a description of the power cepstrum, it 



4 

then discusses the complex cepstrum and the process of incident wave 
recovery. The chapter closes with limitations anticipated when 
applying cepstral analysis to ocean waves. 

Chapters Four and Five discuss the test models used in 
applying cepstral analysis to simulated water waves and the results 
of this application. The first model is a finite summation of sinusoid 
signals containing random phase components keeping both delay time 
and reflection coefficient constant for all frequencies. The second 
model is the same summation of sinusoid signals as the first but the 
delay time of the reflected wave is allowed to vary with frequency. 

The third model is a summation of sinusoid signals allowing the 
reflection coefficient as well as the delay time of the reflected 
wave to vary with frequency. 

Chapter Six is a brief look at some real world ocean data 
using cepstral analysis. A hint is given of problems remaining to 
be solved before this technique can be easily applied in the field. 

Chapter Seven lists the conclusions from this work 
including recommendations for further study. Overall, the process 
shows promise and with more research could become a valuable tool 
in the decomposition of combined incident/reflected sea states. 
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CHAPTER II . 

BACKGROUND INFORMATION 

The work throughout this study is based on use of the 
discrete Fourier Transform, specifically the FFT or Fast Fourier 
Transform. Reference will be made for illustrative purposes to 
other transform methods but all analysis was actually performed using 
the FFT, For a general review of the FFT and what it means and can 
do, the reader is referred to Bergland T s article "A Guided Tour of 
the Fast Fourier Transform" {5}, 

Although for simplicity some examples may be written in 
terminology more akin to continuous transform analysis, all work in 
this study is specifically digital time series analysis and is 
sometimes only possible in a sampled time series sense. Since 
Digital Signal Processing or Digital Time Series Analysis is an 
independent field of study usually couched in the realm of Electrical 
Engineering, strict adherence to digital analysis terminology would 
possibly thwart the main thrust of this paper which is to explain the 
cepstral analysis process in terms an engineer or scientist who is 
not thoroughly familiar with time series analysis may still understand. 
For a detailed description of the Discrete Fourier Transform, the 
reader is referred to Chapter 9 of Stanley {6}, Chapters 3 and 6 of 
Oppenheim and Schafer {7}, Chapters 6-9 of Brigham {8}, and/or 
Chapter 5 of Childers and Durling {9}. 

All computer modeling and testing was performed using CSAP 
(Cepstral and Spectral Analysis Program) which was developed by the 
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author from TSAPC, an instructional program developed by Christensen 
at the University. of Washington used in teaching spectral analysis. 
Some subroutines for CSAP were obtained from NSAP (Numerical Spectral 
Analysis Program) also developed by Christensen and were modified by 
the author as appropriate to meet the needs of this particular 
research. The FFT subroutine used was originally called NLOGN 
developed by Thompson and is described by Paniker {10}. Spectral 
computation is performed using this FFT subroutine and decimation by 
frequency. The cepstral computations were performed using a sub- 
routine developed by the author based heavily on work presented by 
Skinner {11}. 

Host signal processing people address cepstral analysis in 
the literature in terms of the Z-transform. In this study, use is 
made of the Fourier Transform which is exactly the same operation as 
computing the Z-transform on the unit circle. In brief, the Z- 
transform is merely a discrete version of the Laplace transform. 

The Laplace transform (one-sided) can be expressed by 
X(s) = /“x(t)e~ st dt 

whose inverse is 

x(t) = -i— . / X(s)e st ds 
2iri c 

It can be . seen that if s - ia> where to = 2?rf the above relationships 
become the familiar continuous Fourier transform pair. 

The Z-transform is 

X(z) = I x(n)z" n 
n=o 

sAt 



where z = e' 
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If s = im, then the Z-transf orm* is computed on the unit 

circle which yields the discrete Fourier transform, 

X(f) = Z x(n)e“ la)nAt 
n=o 

or more correctly as the finite sum 

TSE-1 -i27rhn 

X(h) = ^ x(n)e N 
n=o 

(where N = total number of sampled points and h represents digitized 
frequency in the same way that n represents digitized time) 
whose inverse is 



1 Ifel 



i2rhn 



x(n) = - l X(h)e N 
N n=o 

In this way, data points representing digital or discrete time series 
operated on by the Fourier transform become data points representing 
a digital frequency series, and the digital frequency series operated 
on by the inverse Fourier transform becomes a digital time series. 

For our purposes this relationship between the Z-transform 
and the Fourier transform is sufficient* The Z-transform is 
extensively used in the literature as it is a more general 
mathematical principle that includes the Fourier transform as only 
a specific part of its domain* For a detailed background on the 
Z-transform the reader is referred to Chapter 4 of Stanley [6}, 
Chapter 2 of Oppenheim and Schafer {7}, Chapter 2 of Childers and 
Durling {9}, and/or Chapter 2 of Gold and Rader {12}. 
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CHAPTER III . 

DESCRIPTION OF CEPSTRAL ANALYSIS 



A. Power Cepstrum and Determination of Delay Times. 

In 1963, Bogert, Healy, and Tukey published their classic 
article {13} on power cepstral analysis which describes the cepstrum 
(pronounced K^p-striim, which violates a few rules of English, but 
sounds better than Sep^-strum) as a tool for determining delay times 
between reflected wave components and incident wave components 
when recording both simultaneously as a combined wave field. The 
simplest model for the combined time series, y(t), containing only 
one reflected wave component or echo is 

y(t) = x(t)+a-x(t-T) 

where x(t) represents the original or incident wave and t represents 
the delay time between the beginning of the incident wave to the 
beginning of the reflected wave. The a term in the above expression 
is the reflection coefficient, which is the ratio of the reflected 
wave height to the incident wave height. 

The Fourier transform of the reflected echo is obtained 
from the Fourier transform of the original wave time series by 
multiplying by ae”^ 7r ^^ T where f is in frequency dimensions (Hz). 
Therefore the Fourier transform of the combined wave time series is 

Y (f ) = X(f){l+ae -27rifT } 

Since the magnitude squared of the Fourier transform represents the 
power spectrum, we have 
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| Y(f ) | 2 = | X (f ) | 2 {l+a 2 +2acos27rfx) 

which represents the power spectrum of the combined wave time series. 



Taking the log of this power spectrum yields 

log|Y(f)| 2 = log|x(f ) | 2 +log{l+a 2 +2acos27rf x} 

log|Y(f)| 2 = log(l+a 2 )|x(f)| 2 +log(14^“ COs2irfx) 

l+a 2 

Using the log series expansion: 



v 2 3 

lo g (1+x) = X^2L-+i_- 
2 3 



,|x|<l 



we can expand the second term into a convergent series thus: 

log(l+— . cos27rfx) = E C z -Q , ( — cos2Trfx)^ 
l+a 2 k=l k l+a 

The log of the power spectrum of the combined wave time 
series then contains a cosinusoidal ripple associated with the 
reflected wave component. The "frequency" of this ripple can be 
found by again calculating a power spectrum. This second power . 
spectrum of the log of the original power spectrum is called the power 
cepstrum. In the transformation process, a change is made from the 
time domain (sec) to a frequency domain (Hz) and then to another 
temporal domain called the quefrency domain (sec) which could be 
thought of as ripples/Hz. 

If we consider a<<l then we can simplify such that 
log|Y(f)| 2 = log(l+a 2 ) |x(f ) [ 2+ ? cos2Trfx) k 



k=l 



1+or 



becomes 



log|Y(f) | 2 = log|x(f ) | 2 +2acos2Trfx 
where the second term in this expression is represented by the 
first term in the convergent series. The power spectrum of this 
expression then yields the power cepstrum 
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ypcOO = Xp C (t)+a6.(t-T) 

such that the power cepstrum of the combined wave series is the 
power cepstrum of the incident wave series plus a unit impulse 
function or spike at the delay time x. In the more general case 
when a is not «1, additional impulses occur at multiples of x. 

If the delay time is sufficiently long, these spikes can 
be observed separated from the incident wave power cepstrum. Thus 
this procedure can be followed for determining the delay time, x, 
between incident and reflected wave. The next step is to try and 
filter out the reflected information found at x and multiples of x 
to attempt recovery of the incident wave. 

B. Complex Cepstrum and Wave Recovery. 

A different approach to cepstral analysis is seen with the 
investigation of Oppenheim and Schafer into homomorphic deconvolution 
and the complex cepstrum {14, 7 (chap. 10)}. Using the single 
reflection model again: 

y (t) = x(t)+ax(t-r) 

The Fourier transform of this is 

Y (f ) = X(f)(l+ae~ 2lTifT ) 

Noting that Y(f) and X(f) are usually complex quantities, 
for this analysis, the log is taken directly of the Fourier 
transform of the time series and not of the power spectrum as was 
done before. We then have 

log{Y (f ) } = log{X (f ) }+log{l+ae“ 2rrifT } 
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Using the log series expansion expressed previously the 
second term becomes 

log(l+ae-2TrifT) = ae ~ i2ir£x -2l 0 -i2(2TrfT) 4Xt 3 0 -i3(2Trfx)_ 

2 3 

The complete log expression then undergoes inverse Fourier 

transformation to yield the complex cepstrum for the combined wave as 

?(t) = X(t)+a6(t-T)-|^5(t-2T)+|is(t-3T)- 

Note that the complex cepstrum for the combined wave time series, y(t), 

consists of the complex cepstrum of the incident wave time series, 

&(t), plus impulse functions or spikes at multiples of x which decay 

exponentially with a and also alternate signs. In all of the above 

a is considered less than 1.0 which is the normal characteristic of 

the reflection coefficient In the ocean. If a is somehow greater 

than 1.0, a similar examination reveals that the complex cepstrum 

becomes the inverse Fourier transform of the complex log of the 

echo, plus impulse functions or spikes which decay exponentially 

with a and alternate signs in negative time only {15}. Kemerait and 

Childers {15} also go through a similar process for the double 

reflection case which yields this complex cepstrum 

£(t) = ^(t)+a 1 (S(t-T 1 )+a 2 6(t-T 2 )-^l_jS(t-2T 1 )-a 1 a 2 iS{t-(T 1 +-r 9 )} 

-fL$(t-2T 2 ) 

2 

where a^, a 2 and x ^ , x 2 represent the different reflection co- 
efficients and delay times for the two different reflections or 
echoes from the same incident wave. Here ja 1 e”^ifT 1+a ^ e -27rifx 2 | 
must be less than 1.0. This could be expanded further to higher 
numbers of multiple echoes from the same incident wave. In the 
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ocean this would be equivalent to recording reflected wave data from 
say two or more neighboring breakwaters simultaneously combined 
with incident wave data. The present effort is limited to the single 
reflection case only. The incident ocean wave can also be thought of 
as consisting of a su mm ation of many different independent waves each 
with its own single reflection but more will be said on this later. 



transform of a real time sequence, some important considerations are 
worthy of mention. First, since with ocean or water wave data, only 
real information is recorded, the Fourier transform of this real 
information yields complex numbered information. These complex 
numbers are manipulated to determine the phase and frequency ampli- 
tude associated with the sampled wave series, which becomes a sampled 
frequency series after transformation. This amplitude is found by 
computing the resultant of each complex vector and the phase, <j>, 
is found by taking the inverse tangent of the imaginary part divided 
by the real part. In the x - iy plane this looks like 



In taking this complex logarithm of a forward Fourier 



^ iy 




r aal part 



imag. 

part 



> x 



Amplitude, |X(f) | 





Phase, <j> 
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The log of these amplitudes (remembering that in discrete 
analysis each data point is acted upon individually) becomes the 
component in the complex cepstrum computation discussed above. For 
wave recovery, however, the phase information cannot be discarded. 
Additionally, since the phase thus computed from complex number 
algorithims is discontinuously varying between -it and -Hr , it must 
be manipulated prior to taking the inverse transform as will be 
discussed later. 

The cepstrum is a real valued quefrency (time) series. 

As such it has no imaginary part. To achieve this, the real part of 
the frequency series undergoing inverse transformation must be con- 
tinuous and even and the imaginary part must be continuous and odd. 
The amplitude of the frequency series computed as shewn above as well 
as the log amplitude will always yield a continuous and even function 
when the original time series is real. The phase, however, must be 
made continuous by "unwrapping” and its linear trend must be removed 
to avoid undesired shifting. When this happens the phase then 
becomes continuous and odd. The log amplitude becomes the real part 
and the unwrapped phase with linear trend removed becomes the 
imaginary part of a frequency series which undergoes inverse Fourier 
transformation to yield the real valued complex cepstrum. 

To explain the phase unwrapping procedure further, Figures 
1 & 2 show a sample raw phase and sample unwrapped phase. The 
unwrapped phase linear trend is removed by adjusting data points so 
that they fall above and below a zero-mean line the proportional dis- 
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tance they occur from a line drawn through the first and last points 
of the unwrapped phase plot. The heavy dotted line on Figure 2 
indicates the line from which these distances from a zero-mean are 
computed. More examples of phase unwrapping and linear trend removal 
will be presented later. 

Skinner {11} has shown that the power cepstrum can be 
computed from the complex cepstrum quite easily as 

y PC ( n ) = (?(n)+?(-n)) 2 

That is, the power cepstrum is merely 4 times the square of the even 
part of the complex cepstrum. He also defines a phase cepstrum that is 
related to the phase in the same way that the power cepstrum is related 
to the log magnitude. The phase cepstrum is four times the square of 
the odd part of the complex cepstrum and looks like 

y_, _(n) = (?(n)-?(-n)) 2 
PhC 

The phase cepstrum is presented only for curiosity as the power 
cepstrum is still a better indicator of the delay time, t. The 
phase cepstrum does, however, give some indication as to whether 
the phase unwrapping procedure was successful. That is, if the 
phase cepstrum is very weak or non existant, the phase unwrapping 
may not have been complete. This could happen if the phase differed 
by more than 2tt between any two sample points. 

The procedure for recovery of the incident wave first 
involves finding the delay time, t, by either the complex or power 
cepstrum. Then the impulse functions or spikes which contain the 
reflected wave information found in the complex cepstrum at t and 
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multiples of t are filtered out and smoothed over by averaging 
neighboring points on either side of the spikes. The resulting 
filtered complex cepstrum then contains information pertaining solely 
to the incident wave. By reversing the process for computing the 
complex cepstrum the original incident wave time series is obtained. 
The algorithm for this process as obtained from Skinner {16} is 
described by Figure 3. 

To demonstrate the wave recovery process, two models are 
taken from Skinner {11} and results plotted using a Calcomp Plotter. 
The first model is of a 256-point time series consisting of two 
64-point long exponential signals with the second one reduced in 
amplitude by one half and delayed 30 points. Using a sampling rate 
of .3333 sec, this means the combined signal is a decaying 
exponential pulse 21.33 sec long combined with its echo reduced in 
half and delayed 10 sec. Mathematically it looks like: 
y(t) = x(t)+.5x(t-x) where x(t) = te -t ; t is actually the sample 
number n multiplied by the sampling rate At and x=30At=30x. 3333 or 
x=10 sec. 

Figure 4 shows the raw spectrum of this combined signal, 
and Figure 5 shows the log amplitude. Figures 6-8 show the raw 
phase, the unwrapped phase, and the phase with linear trend removed. 
Figure 9 shows the real -valued complex cepstrum. Figures 10 and 11 
show the power and phase cepstrums as computed from the complex 
cepstrum. Note the large sharp impulse spikes at multiples of x(10 
sec) in all three cepstral plots and the alternating positive/ 
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negative impulse peaks in the complex cepstrum. Figure 12 shows the 
true incident wave in the middle with the recovered incident wave 
above and the combined wave below. Recovery with this example is 
almost perfect with a Mean Squared Error (MSE) of only 4.44 x 10"^ 
between true incident and recovered incident waves. 

The other signal observed by Skinner {11} was a 256-point 
combined time series containing two FM signals pulsed over 64 points 
of which one was reduced by H in amplitude and delayed by 30 points 
before being added to the incident wave signal. Mathematically 
this combined signal looks like y(t) = x(t)+.5x(t-T) where 
x (t ) - sin((.4+.llt)t). Using a At of .3333 sec, the incident 
FM signal pulse is 21.33 sec long, and T=10sec. Figures 13 through 
21 show the same sequence as for the exponential signal but the 
recovered wave similarity to the incident wave in Figure 21 is much 
more impressive here as the combined signal went through quite a 
bit of unravelling. Note the small aberrations on the recovered 
signal (top time series of Figure 21) where not quite all of the 
reflected wave component was taken out probably due to a bit of 
leakage, etc. associated with the FFT. Mean Square Error (MSE) 
between the recovered incident and true incident waves in this 
case was computed to be 1.91 x 10"^, again indicating almost perfect 
recovery. 

Both of these examples were run under basically noise 
free conditions. Much of the work to date in cepstral analysis 
has been associated with effects of noise added to the signal. 
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Since in the ocean environment wave measurements can be made without 
noise other than the usually negligable amount produced electronically 
during sampling and recording, these effects will not be examined here. 
Generally any noise added to the signal makes recovery that much more 
difficult. 

C. Anticipated Limitations Related to Application with Ocean Waves 

A basic assumption in our original model used to define ceps- 
tral analysis is that the delay time, x, and the reflection coefficient, 
a, are constant and not functions of frequency. With water waves this 
is not a correct assumption. Both x and a may vary with frequency even 
in areas where the water depth is constant. When water depth varies 
and particularly if the waves under study are moving into shallower wa- 
ters, the problem is even more complex due to wave celerity varying with 
depth by different amounts for different wave frequencies. It is hoped 
that the amount of differing x T s and a T s with different frequencies will 
be small enough that this will not cause an insurmountable problem. In 
this regard results should probably be better for single peak and nar- 
row banded spectrums than for double or multi-peaked spec trams with 
broader frequency components as the different individual x r s will hope- 
fully join together to form one main x value. The varying a could pos- 
sibly cause the amplitudes of the impulse spikes in the cepstrum at 
multiples of x to change differently than expected from theory. Both 

£ 

Although the literature shows about equal preference between the use 
of spectra or spectrums as the plural form of spectrum, the latter form 
is used throughout this work. 
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of these effects may require a wider filtering and smoothing of the re- 
flected information spikes in the complex cepstrum to completely remove 
all reflected wave information and allow recovery of the incident wave* 

Another problem that may be encountered using the ocean wave 
is the magnitude of a. Kemerait and Childers {15} found problems try- 
ing to recover the incident wave as the reflection coefficient exceeds 
0.8 for various computer generated test signals. For the multiple re- 
flection case they found that if the reflected signal amplitudes are 
less than the basic wave, recovery is possible. In the general wave 
case, if the reflected signal amplitudes are greater than the incident 
wave, but not equal, the largest reflected signal can be recovered. If 
two or more reflected signals have larger amplitude than that of the 
incident wave, recovery is not possible. Since ocean structures fre- 
quently may have reflection coefficients approaching 1.0, this real 
limitation will be examined for at least the single reflection case. 

Cepstral analysis to date has been involved with pulsed and 
burst-type signals of relatively short duration similar to the recent- 
ly discussed exponential and FM signal examples. A problem may arise 
from using combined incident/reflected wave data collected over the 
total sample length. Since the time series is in fact digitally sam- 
pled it is felt this should not be an insurmountable problem. It may 
mean zeroes will have to be added to the end of the sampled data 
which increases the frequency sampling rate to get results similar to 
that of the pulsed signal. More will be said about this in Chapter V. 



19 



CHAPTER IV 



DESCRIPTION OF TEST MODELS 



A, Summation of Finite Number of Sinusoids with Random Phase 



The first test model for ocean simulated waves maintains 



both a, the reflection coefficient, and t, the delay time, as 
constants non- variant with respect to frequency. The incident wave 
model used in this test is that derived by Borgman {17} which simulates 
an ocean waveform based on an assumed spectral curve. This model can 



having random phase shift angles. Simplified mathematically it looks 



where x(t) represents the simulated ocean wave surface profile. As M, 
which represents the number of sinusoid components, gets larger, the 
more realistic the simulated wave becomes. The tradeoff for using 
larger and larger values of M is increased computer cost for "building" 
the simulated wave. For the testing in this work, an M value of 40 is 
used for a 256-point time series. The amplitude of each component is 
found by dividing the total assumed energy represented as the var- 
iance, a 2 , of the assumed spectrum (where a 2 is the area under that 
spectral curve) by the number of sinusoids being summed and taking the 
square root of that result (remembering that the square root of 
energy is directly proportional to wave height). Therefore, the total 



be described as the summation of a finite number of sinusoid components 



like 




energy in the simulated wave will equal the energy represented by the 
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initially assumed spectral curve. When the assumed spectral curve is 
divided into equal area increments, f m values are chosen as the mid- 
frequency of each incremental frequency band. At is the time increment 
between samples, and <f> m is a uniformly distributed random variable 
between zero and 2tt which represents the random phase shift or 
starting point of each sinusoid component being summed. 

For example, with an M value of 5 the assumed spectrum would 
be divided as shown below to obtain the equal energy value for comput- 
ing amplitudes and for determining f values or the mid-point fre- 



that proposed by Pierson and Moskovitz {18} and modified for limited 
fetches by Mitsuyasu {19}, Barnett {20} and Sylvester and 
Vorgvisessomjai {21}. The general form of the P-M spectrum is 



quencies . 




frequency 



(not to 
scale) 



f f f f f 
1 2 3 5 



The assumed spectrum used in creating this wave form is 



S(f) - 




where a (not to be confused with the reflection coefficient) 



and 3 are assumed to be constant for a fully arisen sea state. U 



is the wind velocity, f the frequency, and g the gravitational 
constant. However, as derived by Christensen and described by 
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Morden {2}, for limited fetches a and 3 become 

a . SMiQVfJ; „ 3.49(2a) l> S(f I ,)f p 5 

P g 2 

s 

where f is the frequency of the spectral peak and S(fp) is the 
amplitude of the spectral peak. The variance a 2 is then found from 



a 2 = 2/”s(f)df = aU . 

° 4gg 2 



A model wave then can be created based on the variables U, 
fp, S(fp), and a 2 . Assuming U and f and either a 2 or S(fp) will 
completely define the P-M spectrum. 

Therefore, the simulated incident ocean wave is computed 
by first assuming values of S(fp) and f p to compute the P-M spectral 
model and a value for a 2 . This spectral curve is divided into equal 
areas and mid-point frequencies for each equal area increment are 
computed. Using a suitable random number generator to compute random 
phase angles, each sinusoid is computed and summed to the previously 
computed sinusoid (s) as described previously. For simplicity, this 
simulated incident ocean wave will be called the P-M generated wave. 

Taking the P-M generated wave and scaling it by a 
reflection coefficient a produces the simulated reflected wave. 
Delaying this reflected wave by x and adding the two together yields 
the initial P-M generated combined wave test model. Note x is 
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constant for the whole wave train in this initial model. In Chapter V, 
further description of model manipulation is described related to test- 
ing procedure. 



B. Varying Delay Time with Respect to Frequency 



In water waves, the delay time t varies with respect to the 

wave celerity which varies as a function of wave length which in turn 

varies with frequency, so it follows that x does in fact vary with 

frequency. For linear deep water wave theory this variation can be 

xf 



expressed in fps units by x 



where f is the frequency of the 



2.56 

wave component and x is the distance from the wave measurement device 



to the reflecting surface. This relationship is derived from 

distance traveled 

celerity = : Z I = 

time of travel 

where T is the wave period and g is the gravitational constant. 




T for deep-water wave theory 



2x 1 

Therefore celerity 3 — - 5.12 — since distance x is traveled twice 

x f 

over the delay time, x, so that x = 

5 • 12 ^ 

Using the P-M generated wave model previously described, the 
incident wave is generated in the same manner as before. However, the 
reflected wave components are delayed by their own x values for their 
particular frequency before being summed together, scaled by a and 
combined with the incident wave. This is then a model where each sep- 
arate reflected frequency component has its own delay time x contri- 
buting to the overall resultant x of the reflected wave time series. 



Results from using this model are discussed in the next chapter. 
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C. Varying Reflection Coefficient with Respect to Frequency, 

To carry the P-M generated wave model one step closer to 

the real world reflected water wave, an attempt is made at varying 

the reflection coefficient as a function of frequency. This first- 

try variation is performed by merely making a a linear function vary- 

. 2 

ing from .4 to .6, .5 to .7, or .7 to .9 incremented by for each 
sinusoidal component of the reflected wave where M is the total number 
of components. The results from using this test model are also dis- 
cussed in Chapter V. 

In reality, the reflection coefficient, a, probably varies 
as some non-linear function of frequency such that ct would be higher 
for higher frequency waves than for lower frequency waves depending on 
the structure. That is, the lower frequency wave components are more 
likely to be absorbed or transmitted by most structures than are the 
higher frequency wave components making up the ocean wave. However, 
the simplified model described above and in Chapter V is sufficient 
for an initial test of effects from ct varying as a function of 
frequency. 



24 



CHAPTER V 

APPLICATION OF CEPSTRAL ANALYSIS TO TEST MODELS 

As mentioned in Chapter III, cepstral analysis techniques 
are usually applied to "pulse" type signals where no signal informa- 
tion is contained between the incident wave and the reflected wave. 
Figure 22 shows a sample of incident and reflected pulse signals which 
would be added together to form such a combined incident /reflected sig- 
nal* Note the lack of information over time, t, at the end of the 
incident wave and over time, x, at the beginning of the reflected wave* 
Both incident and reflected waves are of the same length. The expon- 
ential and FM signals discussed in Chapter III are examples of pulse 
waves. 

When the incident wave extends the full length of the time 
series, the incident and reflected waves appear as shown in Figure 23. 
This combined wave is called a "pulse-shifted" wave as no infor- 
mation exists over time, x, in the beginning of the reflected wave. 

Note the reflected wave is now not as long as the incident wave. This 
creates a situation known as "echo truncation" where the incident wave 
does not terminate x seconds before the end of the reflected wave as 
in the pulse signal described above. 

For a more realistic look at ocean wave data sampling, a 
"forward-shift" combined wave is created from incident and reflected 
waves following the pattern shown in Figure 24. The name "forward- 
shift" comes from the fact that a new incident wave is made from the 
previously created incident wave by extending the wave train and 
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figuratively shifting it forward. The important item here is that 
signal information exists from both the incident and reflected waves 
throughout the sampling period. 

To test the first model described in Chapter IV where neither 
a nor t are allowed to vary with frequency the combined wave is formed 
first using the "pulse-shift” combination and then using the "forward- 
shift" combination. In both cases, the incident and reflected waves 
are created and stored and each data point in the reflected wave is 
multiplied by a reflection coefficient before the two waves are added 
together to form the combined wave. 

The P-M Spectral model used to create the P-M Generated Wave 
is shown as Figure 25. This spectrum is selected as most of the spec- 
tral energy is concentrated around a narrow frequency band. With later 
trials this should help produce one x value for the total combined wave 
when x is allowed to vary with frequency. 

Using a rectangular window (no weighting of input data) as 
recommended by Skinner {11} for the pulse-shifted wave, results are 
somewhat discouraging. For a 256 data point time series. At-. 3333 sec, 
x=10 sec (30 data points) and a=.8, the Mean Squared Error (MSE) between 
the recovered incident wave and the true incident wave is .122 with the 
MSE between the true incident and combined wave equal to .285. The MSE 
(equal to the mean of the difference squared comparing parallel data 
points of two "side-by-side" time series) between recovered and true 
incident waves is an indicator of recovery results and the MSE between 
true incident and combined waves is an indicator of what the recovery 
process starts with. Comparing these two indicators, the MSE relation- 
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ships show that some recovery did occur but that the recovered incident 
wave and the true incident wave are still not very much alike. 

Oppenheim and Shafer {7} state that using an exponential win- 
dow (exponential weighting of input data) enhances the identification 
of delay time, t, and improves the overall cepstral analysis recovery 
process. Skinner {11} found this not to be the case for the exponen- 
tial and FM signals discussed in Chapter III. Using an exponential 
window on the pulse-shifted combined wave where each data point (n=l to 
n=N) is "weighted" or multiplied by (.99) n ~^ before the cepstral decom- 
position process and dividing the resulting recovered wave by (.99) n 
to "unweight" the output lowers the MSE between original and recovered 
incident waves to .044. This indicates quite satisfactory recovery re- 
sults as can be observed from Figure 34. 

It should be pointed out here that although "one-sided" plots 
(representing only the positive half of the spectrum) of the raw spec- 
trum and log amplitude are shown for the exponential and FM signals, 
all subsequent plotting of raw spectrum and log amplitude before the 
cepstrum is computed are shown as "two-sided" plots (noting that the 
negative half of the spectrum is reproduced to the right of the posi- 
tive half since transformation was performed via the DFT) . This is 
done to show that in fact the log amplitude is indeed continuous and 
even and the unwrapped phase with linear trend removed is continuous 
and odd which together produce a real valued cepstrum. Smoothed spec- 
tral curves to be seen later will still be shown as one-sided plots. 

To discuss the positive results from exponential weighting 
further, it is seen that it helps reduce echo truncation error by plac- 
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ing greater emphasis on the beginning of the input combined wave time 
series than on the end. Different exponential "base" values from .900 
to .999 were tried for different input conditions and a base value of 
.99 always seemed to produce the best results. 

Figures 26 through 34 show the plots of the time series 
manipulation through. the spectral and cepstral processes at the same 
stages as listed for the exponential and FM signals. An important ob- 
servation when seeing these plots is that this more irregular wave form 
produces a more irregular cepstrum but the reflected information im- 
pulse spikes are sti],l very obvious in all three cepstrums. When the 
spikes located at multiples of t are filtered out of the complex ceps- 
trum, wave recovery is quite good as is seen in Figure 34 when compar- 
ing the original incident time series located in the middle of this 
figure to the recovered incident time series located at the top . 

A better indicator of how good the results are using this 
model is seen in Figures 35 through 38. These figures represent the 
smoothed spectrums (every 8 points averaged together for these parti- 
cular plots) of the combined wave, the original incident wave, the re- 
covered incident wave, and the reflected wave, respectively. The spec- 
trums of the recovered and original incident waves are almost identical 
indicating good wave recovery. The reflected wave spectrum is obtained 
from the recovered reflected wave found by subtracting the recovered 
incident wave data points from the combined wave data points. Using 
the formula which relates the variance, or area under the 

spectrum of the reflected wave to the variance or area under the spec- 
trum of the incident wave times the reflection coefficient squared, an 
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attempt is made at determining the reflection coefficient used to cre- 
ate the combined wave. The resulting computation produces a reflection 
coefficient equal to .785, which is very close to the actual .8 value 
used. 

The results of this first test of a P-M generated wave shape 
therefore are very encouraging. The next tests using the "forward- 
shift” model, however, are not so gratifying. Input parameters are the 
same as with the pulse-shift model. Figures 39 through 47 show the 
same plot sequence as before. The impulse spikes are not as readily 
noticable in the three cepstrums now and recovery attempts as seen from 
the time series comparisons in Figure 47 between the original incident 
wave (middle) and the recovered incident wave (top) show results poorer 
than with the pulse-shift test. The smoothed (averaged spectral value 
of every 10 data points) spectral plots of results shown in Figures 48 
through 51 also show recovery as less than satisfactory. The recovered 
wave spectrum shows its energy level to be somewhere between that of 
the true incident and combined waves instead of equal to that of the 
true incident wave. The reflected wave spectrum is obviously too small 
such that the computed reflection coefficient is much lower than the 
value of .8 used in creating the combined wave. An encouraging obser- 
vation in this, however, is that some energy is in fact removed from 
the combined wave so that the process at least attempted to perform 
some decomposition. 

Following Skinner’s {11} recommendation, a first step or try 
at improving results is to add N zeroes to the end of the input data 
(adding 256 zeroes for 256 data points) before performing the spectral 
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and cepstral computations. This effectively doubles the sampling rate 
in the frequency domain and hence the quefrency domain without changing 
the time series information whatsoever. When zeroes are added to the 
end of the pulse-shifted combined wave from before, no improvement is 
noted, probably due to the fact that results are already very good. 
Adding N zeroes to the end of the input data of the forward-shifted 
combined wave before performing decomposition yields the smoothed re- 
covered incident and reflected wave spectral plots shown in Figures 
52 and 53. It is seen that this technique helps produce a recovered 
incident wave spectrum that is slightly closer to the true incident 
wave spectrum, but results are still unsatisfactory as evidenced by the 
again small reflected wave spectrum and the obvious resulting incorrect 
calculation of the reflection coefficient. 

Manipulation of data sampling methods is also employed in 
searching for better results. These methods include trying different 
time sampling rates, trying different t values and increasing the total 
sample length to 1024 data points from the previous length of 256 data 
points. All attempts at data sampling manipulation seem to be ineffec- 
tual for improving results. When the reflection coefficient is reduced 
to .5, results improve only very slightly. Attempts at improving re- 
covery using a rectangular window and an exponential window with dif- 
fering base values above and below .99 show best results are still 
achieved with exponential weighting and un-weighting using the (.99) n “^- 
factor described previously. 

Therefore, a problem exists using the P-M Generated wave 
forward-shift model even when t and a are constant and independent of 
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frequency. This problem appears to be centered around the fact that 
information is now contained in the beginning of the reflected wave. 
This information looks more like noise when comparing the incident and 
reflected waves’ time series because the actual image reflection from 
the incident wave doesn't start until t seconds later. 

When the second model described in Chapter IV is employed, 
where x Is allowed to vary with frequency for each sinusoid making up 
the P-M generated combined wave, even less incident wave recovery is 
possible. The different attempts at improving results follow the 
same sequence as those tried with the forward shift wave. Figure 54 
shows a typical complex cepstrum for this model where At=.3333, x=10 
sec and a=.5. It is readily obvious that identifying a x, even when 
the input x is known, is very difficult in this situation. Filtering 
the impulse spikes which should be present at multiples of x yields 
the filtered complex cepstrum shown in Figure 55. Figures 56 and 57 
demonstrate the same technique showing a complex cepstrum and a 
filtered complex cepstrum for At=.7333. Again it is quite difficult 
to identify ax. 

Since it appears initially that the x peak in Figure 54 may 
be wider than just one point, five different filter variations are at- 
tempted filtering out more than one point. Some examples are shown in 
Figures 58 and 59 of the resulting filtered complex cepstrums. None 
of these different filter methods help and in fact all of them make 
the MSE even greater between the original and recovered incident waves . 
Some attempts even made the MSE between actual incident and recovered 
incident waves higher than between the incident and combined waves . 
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Therefore, it appears the impulse spike definition for location of x is 
more accurate than a wider peak definition of t even in the model where 
x is allowed to vary with frequency. 

When employing the third model described in Chapter IV, where 
reflection coefficient is also allowed to vary with frequency, very 
little difference is observed in the appearance of the complex cepstrum 
using the same input parameters as for the constant reflection coef- 
ficient model above. Wave recovery attempts for this model are not sa- 
tisfactory. Again, no improvement is obtained by using different At T s v 
different x f s, or different reflection coefficient ranges, including a 
varying from .7 to .9, .5 to .7, and .4 to .6. Example complex ceps- 
trums for this model are shown in Figures 60 and 61. The difficulty in 
identifying a value for t is also readily apparent for this model. 

The problems with the latter two models cannot be readily 
identified until satisfactory results can first be obtained from the 
forward shift model where a and t are non-variant with frequency. 

This requires some way of satisfactorily dealing with the "noise" or 
continuous wave section existing in front of the reflected wave con- 



tribution. 
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CHAPTER VI 

APPLICATION OF CEPSTRAL ANALYSIS TO ACTUAL OCEAN DATA 

Some time series data records recorded near a breakwater in 
Sitka, Alaska, in 1973 are run through this algorithm to determine if a 
value for t can be identified. Figure 62 shows the breakwater plan and 
the location of instruments on this breakwater. The data points util- 
ized are taken from the combined wave recorded at the Northernmost 
gage when wave action was predominately perpendicular to the short leg 
of the breakwater. This site has the advantage that any recovered in- 
cident wave obtained through the decomposition process can be compared 
for accuracy with the true incident wave recorded from the wave staff 
off the corner of the breakwater. 

These data records had a major disadvantage, however, as 
the breakwater is only about two feet deep allowing some of the low 
frequency waves to transmit through the breakwater causing the reflec- 
ted wave to be different in shape and frequency composition from the 
incident wave. Also, the frequency spectrum for this location is 
’'double-peaked,” that is, a low frequency swell condition exists as 
well as the waveform produced by local winds over a local fetch. 

Looking at cepstrums generated from this data is informative, 
however, in that it reveals some of the other problems associated with 
using this process on real world data. Three records are selected and 
tested which are 2048 data points long recorded at a At of .44 sec. 
producing 15 minute long data records. An example of the complex 
cepstrum thus obtained using no weighting of the input data 
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(rectangular window) is shown in Figure 63.. It is interesting 
to observe the unwrapped phase and unwrapped phase with linear trend 
removed for this data which are shown in Figures 64 and 65. The 
similarity between the two raises suspicion that phase unwrapping 
may not have been complete. A remedy for this problem is to add 
zeroes to the end of the original input data. This is tried using 
1024 data points from the same record and adding 1024 zeroes before 
computing the spectrum and complex cepstrum. 

Figures 66, 67, 68, and 69 show the raw spectrum, log 
amplitude, unwrapped phase, and unwrapped phase with linear trend 
removed for this new input. The linearity seen in the unwrapped phase 
plot shows phase unwrapping is indeed now more complete. 

Figures 70 and 71 show the complex cepstrum and power 
cepstrums, respectively. Difficulty in obtaining a value for x from 
these plots is readily apparent. Perhaps some help could be obtained 
by plotting the cepstrum over a longer quefrency axis for better plot 
visualization. Using deep water wave theory, the physical location 
of the wave staff and predominant frequency from the wave spectrum, 
t is estimated to be about 18-20 seconds. This would probably be 
very hard to identify in a 15 minute record even if an impulse spike 
were located at this point. 

The most important reason for failure to identify a t for 
this situation is more probably due to the fact that the reflected 
wave is altered in shape and frequency and is not in fact a true 
replica of the original incident wave. It seems Cepstral Analysis 
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decomposition should work even if the reflected wave shape is 
different from the incident wave shape as long as the frequency 
information is not altered in the reflection process. Since the 
floating breakwater used here does in fact filter certain low 
frequency components out of the reflected wave by allowing them to 
pass through the breakwater, it is not surprising that an impulse 
spike is not observable in either the complex or power cepstrums at 
an identifiable location for x. 

The predominant spike at the beginning of these cepstrums 
is left over from the incident wave complex cepstrum which always 
produces an extremely large spike at zero quefrency in the cepstral 
plots. In fact, all of the cepstral plots seen in this thesis have 
their first 5 data points set to zero amplitude on either side of 
zero quefrency (more properly in the digital sense, points 
n ■ 1, 2, 3, 4, 5, N, N-l, N-2, N-3, and N-4 were set equal to zero) 
to allow easier identification of the smaller amplitude spikes 
containing the reflected wave information. Otherwise, plotting 
would be based on the large spike at zero quefrency and the smaller 
reflected information impulse spikes would be even more difficult 
to identify. 

Even taking this into consideration, an actual x value 
is practically impossible to identify from these cepstrums. The 
power cepstrum shown in Figure 71 offers many possible impulse spikes 
but at times which are totally unrealistic with the ocean wave data 
interacting with this breakwater. Obviously, better data and more 



research are required to demonstrate any sort of positive results 
with cepstral analysis decomposition 



using ocean data. 
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CHAPTER VII 

CONCLUSIONS AND RECOMMENDATIONS 

Cepstral analysis decomposition techniques work extremely 
well with the pulse-shift P-M generated combined wave model and not 
so well with the forward-shift P-M generated combined wave model as 
described in Chapters TV and V. Problems exist with identification 
of delay time as well as with decomposition of the combined wave when 
cepstral analysis techniques are applied to P-M generated wave 
models which vary x and/or the reflection coefficient, a, as a 
function of frequency. More research is required before this 
decomposition technique can be successfully applied to real-world 
ocean waves. It may never work for decomposition of combined 
waves for which the reflected wave frequency information is 
altered by the structure and no longer is the same as the incident 
wave, such as occurs with shorelines and many breakwaters. 

Exponential windowing or weighting does help when using 
these decomposition techniques with P-M generated wave models. 

Adding zeroes to the end of the input data prior to application of 
cepstral analysis decomposition does improve phase unwrapping but 
helps only slightly to none at all with incident wave recovery. 

Of the three cepstrum types examined, the power cepstrum 
is the best .indicator of delay time between incident and reflected 
waves. The phase cepstrum does not seem to offer much information 
other than reinforcement or possible validation of the power 
cepstrum information. 
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It is recommended that this process be tested using real- 
world ocean data for which the reflected wave frequency information 
is identical to that of the incident wave such that all frequency 
components are reflected by the structure. Further experimentation 
should begin using real-world data with a "single-peaked" smoothed 
spectrum before going to data with a "two-peaked" or "multi-peaked" 
spectrum. When using data with more than one spectral peak, 
experimentation could be made into the effect of filtering out one 
or more peak(s) before going through the decomposition process. Of 
particular interest would be results obtained after filtering out 
the lower frequencies from the combined wave that would be 
transmitted through a structure like the floating breakwater 
investigated in Chapter VI. 

Since the complex cepstrum consists of impulse spikes 
at multiples of x which decay exponentially as a function of the 
reflection coefficient, investigation could be made into possible 
determination of reflection coefficients directly from the 
amplitudes of these peaks. This would save the inverse process 
presently necessary by not requiring actual determination of the 
incident wave to obtain the reflected wave before computing a 
reflection coefficient. 

Research should be conducted towards improving results 
with the forward-shift. One possible alternative is to try 
reiterative techniques. Once a value for delay time is determined, 
the recovered wave could possibly be improved by putting it through 
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the decomposition process again. Another possibility to improve 
results is to somehow weight or alter the data points at the 
beginning of the combined wave time series so that the leading 
portion of the reflected wave will not have such a degrading effect 
on results. 

It appears from attempts at adjusting the filtering of 
the impulse spikes from the complex cepstrum that the spikes do 
not join together to form a wider peak when x varies with frequency. 
Perhaps what actually happens is impulse spikes occur separately 
at the different x values produced by different frequency 
components making up the wave. This possibility could be 
investigated further. 

Other windowing or weighting of input data prior to the 
decomposition process has generally not met with success when 
using cepstral analysis techniques. However, this could be 
investigated further for the ocean wave by experimenting with other 
data windows. 

Generally, cepstral analysis decomposition methods do not 
appear at this time to be the final answer in recovery of the 
incident wave from combined incident/reflected ocean wave data. 

First, the x value itself must be made reliably discernable in the 
cepstrum and then actual recovery must be made more efficient. 

However, with further research, it may yet prove to be a valuable tool 
in helping to determine delay times, reflection coefficients, and 
other parameters associated with wave decomposition when determining 
true wave forces acting on and/or from a structure. 



1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 

13 



39 



REFERENCES - 

Ippen, et al., Estuary and Coastline Hydrodynamics, McGraw Hill 
Book Co., New York, 1966, pp. 382-385. 

Morden, D.B. "The Decomposition of Co-Existing Random Incident 
and Reflected Wave Energy," PhD Dissertation, University of 
Washington, 1975. 

Morden, D.B., Richey, E.P. & Christensen, D.R., "Decomposition 
of Co-Existing Random Wave Energy," Proceedings of the Fifteenth 
Coastal Engineering Conference , 1976, pp. 846-865. 

Thornton, E.B., & Calhoun, R.J., "Spectral Resolution of Break- 
water Reflected Waves," ASCE Proceedings, Journal of Waterways, 
Harbors, and Coastal Engineering Division, Vol. 98, November, 
1972, pp. 443-460. 

Bergland, G.D., "A Guided Tour of the Fast Fourier Transform," 
IEEE Spectrum , Vol. 6, No. 7, July 1969, pp. 41-52. 

Stanley, W.D. Digital Signal Processing , Reston Publishing Co., 
Reston, Virginia, 1975. 

Oppenheim, A.V., and Schafer, R.W., Digital Signal Processing, 
Prentice Hall, Englewood Cliffs, New Jersey, 1975. 

Brigham, E.L., The Fast Fourier Transform, Prentice Hall, Inc., 
Englewood Cliffs, New Jersey, 1974. 

Childers, D., and Durling, A., Digital Filtering and Signal 
Processing, West Publishing Co., St. Paul, Minn., 1975. 

Paniker, N.N., "Determination of Directional Spectra of Ocean 
Waves from Gage Arrays, " University of California Hydraulic 
Engineering Laboratory Research Report HEL 1-18, August, 1971. 

Skinner, D., "Real Time Composite Signal Decomposition," PhD 
Dissertation, University of Florida, 1974. 

Gold, B., and Rader, C.M., Digital Processing of Signals, McGraw 
Hill Book Co., New York, 1969. 

Bogert, B. , Healy, M., and Tukey, J., Quefrency Alanysis of 
Time Series for Echoes: Cepstrum, Pseudo-Autocovariance. 

Cross-Cepstrum and Saphe Cracking , M. Rosenblatt, Ed., New York, 
Wiley, 1963, Ch. 15, pp. 209-243. 



40 



14 Oppenheim, A.V., Schafer, R.W., and Stockham, J.G., "Nonlinear 
Filtering of Multiplied and Convolved Signals," Proceedings of 
the IEEE, Vol. 56, No. 8, August, 1968, pp. 1264-1291. 

15 Kemerait, R.C., and Childers, D.G., "Signal Detection and 
Extraction by Cepstrum Techniques," IEEE Transactions on 
Information Theory , Vol. IT-18, No. 6, November, 1972. 

16 Skinner, D.P. and Childers, D.G. "Real-Time Composite Signal 
Decomposition," IEEE Transactions on Acoustics, Speech, and 
Signal Processing, June, 1976, pp. 267-270. 

17 Borgman, L., "Ocean Wave Simulation for Engineering Design", 
Journal of the Waterways and Harbors Division, Proceedings of 
the ASCE , November, 1969, pp. 557-583. 

18 Pierson, W.J. and Moskovitz, L.I., "A Proposed Spectral Form for 
Fully Developed Seas Based on the Similarity Theory of S.A. 
Kitaigorodskii, 11 Journal of Geophysical Research, Vol. 69, 1964, 
pp. 5181-5190. 

19 Mitsuyasu, "One Dimensional Wave Spectra at Limited Fetch," 
Proceedings of the Thirteenth Coastal Engineering Conference , 
ASCE , New York, N.Y., 1972, pp. 351-360. 

20 Barnett, T.P. "Observations of Wind Waves Generation and 

Dissipation in the North Sea: Implications for Offshore 

Industry," Offshore Technology Conference, Paper // 0TC 1516, 
Dallas, Texas, 1972, p. 10. 

Silvester, R. and Vorgvisessomjai, S., "Energy Distribution 
Curves of Developing and Fully Arisen Seas," Journal of 
Hydraulic Research, Vol. 8, No. 4, Delft, Netherlands, 1970, 
pp. 493-521. 



21 



41 



BIBLIOGRAPHY . 

Barber, N.F., Experimental Co rrelo grams and Fourier Transforms, 
Fergamon Press, New York, 1961. 

Barnett, T.P. “Observations of Wind Waves Generation and 
Dissipation in the North Sea: Implications for Offshore Industry,” 

Offshore Technology Conference , Paper #OTC 1516, Dallas, Texas, 1972, 

p. 10. 

Bendat, Julius S. and Piersol, Allan G., Random Data Analysis and 
Measurement Procedures , Wiley-Interscience, a Division of John Wiley 
& Sons, Inc., N.Y., 1971. 

Bergland, G.D., ”A Guided Tour of the Past Fourier Transform,” IEEE 
Spectrum, Vol. 6, No. 7, July 1969, pp. 41-52. 

Bloomfield, Peter, Fourier Analysis of Time Series: An Introduction , 

John Wiley & Sons, N.Y., 1976. 

Bogert, B., Healy, M., and Tukey, J., Quefrency Alanysis of Time 
Series for Echoes: Cepstrum, Pseudo -Autocovariance, Cross -Ceos t rum 

and Saphe Cracking , M. Rosenblatt, Ed., New York, Wiley, 1963, Ch. 15, 
pp. 209-243. 

Borgman, L., "Ocean Wave Simulation for Engineering Design,” Journal 
of the Waterways and Harbors Division, Proceedings of the ASCE , 
November, 1969, pp. 557-583. 

Bracewell, Ron, The Fourier Transform and Its Applications , McGraw 
Hill Book Co., New York, 1965. 

Brigham, E.L., The Fast Fourier Transform , Prentice Hall, Inc., 
Englewood Cliffs, New Jersey,- 1974. 

Brillinger, David R. , Time Series Data Analysis and Theory , Holt, 
Rinehart and Winston, New York, 1975. 

Carnahan, Brice, Luther, H.A., and Wilkes, James 0., Applied 
Numerical Methods, John Wiley & Sons, Inc., New York, 1969. 

Champ eney , D . C . , Fourier Transforms and Their Physical Applications , 
Academic Press, London and N.Y., 1973. 

Childers, D., and Durling, A., Digital Filtering and Signal 
Processing , West Publishing Co., St. Paul, Minn., 1975. 

Churchill, Ruel V., Fourier Series and Boundary Value Problems , 
McGraw-Hill Book Co., 1969. 



42 



Churchill, R.V., Brown, J.W., Verhey, R.F.., Complex Variables and 
Applications , McGraw-Hill Book Co., New York, 1976. 

Cooper, George R., and McGillen, Clare D., Methods of Signal and 
System Analysis , Holt, Rinehart and Winston, Inc., N.Y., 1967. 

Enochson, Loren D., and Otnes, Robert K. , Programming and Analysis 
for Digital Time Series Data, The Shock and Vibration Information 
Center, U.S. Dept, of Defense, 1968. 

Gold, B., and Rader, C.M., Digital Processing of Signals , McGraw- 
Hill Book Co., New York, 1969. 

Goldberg, Richard R., Fourier Transforms , Cambridge at the University 
Press, 1961. 

Hannon, E.J., Multiple Time Series , John Wiley and Sons, Inc., 

New York, 1970. 

Ippen, et al.. Estuary and Coastline Hydrodynamics , McGraw-Hill Book 
Co., New York, 1966„ 

Jenkins, Gwilyn M., and Watts, Donald G., Spectral Analysis and its 
Applications , Holden-Day, San Francisco, 1969. 

Kemerait, R.C., and Childers, D.G., ff Signal Detection and Extraction 
by Cepstrum Techniques, 11 IEEE Transactions on Information Theory , 

Vol . IT-18, No. 6, November, 1972. 

Koopmans, L.H., The Spectral Analysis of Time Series , Academic Press, 
New York, 1974. 

Lubkin, Y.J., Filter Systems and Design: Electrical, Microwave, 

and Digital , Addison-Wesley Bub. Co., Reading, Mass., 1970. 

Mikhail, E.M., Ackermann, F., Observations and Least Squares , IEP- 
a Dun-Donnelley' Publisher , N.Y., 1976. 

Mitsuyasu, M 0ne Dimensional Wave Spectra at Limited Fetch," 
Proceedings of the Thirteenth Coastal Engineering Conference, ASCE , 
New York, N.Y., 1972, pp . 351-360. 

Morden, D.B. "The Decomposition of Co-Existing Random Incident and 
Reflected Wave Energy," PhD Dissertation, University of Washington, 
1975. 

Morden, D.B., Richey, E.P. & Christensen, D.R., "Decomposition of 
Co-Existing Random Wave Energy," Proceedings of the Fifteenth 
Coastal Engineering Conference , 1976, pp . 846-865. 



43 



Oppenheim, A.V., Schafer, R.W., and Stockham, J.G., f, Nonlinear 
Filtering of Multiplied and Convolved Signals," Proceedings of the 
IEEE , Vol. 56, No, 8, August, 1968, pp. 1264-1291. 

Oppenheim, Ed., A.V., Papers on Digital Signal Processing , The M.I.T. 
Press, Cambridge, Mass., 1969. 

Oppenheim, A.V., and Schafer, R.W., Digital Signal Processing , 

Prentice Hall, Englewood Cliffs, New Jersey, 1975. 

Otnes, R.K. and Enochson, Loren, Digital Time Series Analysis , A 
Wiley-Interscience Publication, John Wiley & Sons, New York, 1972. 

Paley, R.E.A.C., Weiner, N., Fourier Transforms in the Complex Domain , 
American Mathematical Society, New York, 1934. 

Paniker, N.N., "Determination of Directional Spectra of Ocean Waves 
from Gage Arrays," University of California Hydraulic Engineering 
Laboratory Research Report HEL 1-18, August, 1971. 

Pierson, W.J. and Moskovitz, L.I., "A Proposed Spectral Form for 
Fully Developed Seas Based on the Similarity Theory of S.A. 
Kitaigorodskii, Journal of Geophysical Research , Vol. 69, 1964, 
pp. 5181-5190. 

Rabiner, Lawrence R., and Rader, Charles M. , editors. Digital Signal 
Processing , IEEE Press, Institute of Electrical and Electronic 
Engineers, Inc., N.Y., 1972. 

Rosenblatt, M. , Ed., Time Series Analysis , Wiley & Sons, N.Y., 1963. 

Silvester, R. and Vorgvisessomjai, S., "Energy Distribution Curves 
of Developing and Fully Arisen Seas," Journal of Hydraulic Research , 
Vol. 8, No. 4, Delft, Netherlands, 1970, pp. 493-521. 

Skinner, D., "Real Time Composite Signal Decomposition," PhD 
Dissertation, University of Florida, 1974. 

Skinner, D.P. and Childers, D.G., "Real-Time Composite Signal 
Decomposition," IEEE Transactions on Acoustics, Speech, and Signal 
Processing, June, 1976, pp. 267-270. 

Stanley, W.D., Digital Signal Processing , Reston Publishing Co., 

Reston, Virginia, 1975. 

Thornton, E.B., & Calhoun, R.J., "Spectral Resolution of Breakwater 
Reflected Waves," ASCE Proceedings, Journal of Waterways, Harbors , 
and Coastal Engineering Division , Vol. 98, November, 1972, pp. 443-460. 



44 



U.S. Army Coastal Engineering Research Center 
Manual, Volumes I, II, III, Superintendent of 
Government Printing Office, Washington, D.C., 



Shore Protection 
Documents, U. S. 
1973. 



Wason, M.T., Parametric Estimation, McGraw-Hill Book Co., 



N.Y., 1970. 



45 




Figure 1. Sample Raw Phase Plot 




Figure 2. Sample Unwrapped Phase Plot 
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Figure 3, Wave Recovery Algorithm 
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Figure 4. Raw Spectrum of Combined Exponential Signal 
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Figure 5. Log Amplitude of Combined Exponential Signal 
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Figure 7. Unwrapped Phase of Combined Exponential Signal 
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Figure 8. Unwrapped Phase with Linear Trend Removed of Combined 
Exponential Signal 
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Figure 9. Complex Cepstrum for Combined Exponential Signal 
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Figure 10. Power Cepstrum for Combined Exponential Signal 
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Figure 11. Phase Cepstrum for Combined Exponential Signal 
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Figure 12. Recovered Incident, True Incident and Combined Wave Time Series for Exponential 
Signal 




Figure 13, 
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Raw Spectrum of Combined FM Signal 
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Figure 14. Log Amplitude of Combined FM Signal 
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Figure 15. Raw Phase of Combined FM Signal 
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Figure 16* Unwrapped Phase of Combined FM Signal 
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Figure 18. Complex Cepstrum for Combined FM Signal 
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Power Cepstrum for Combined FM Signal 
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Figure 19. 
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Figure 20. Phase Cepstrum for Combined FM Signal 
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Figure 21. Recovered Incident, True Incident and Combined Wave Time Series for FM Signal 
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Figure 22. 



Sample Pulse Incident and Reflected Signal 
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Figure 23. Sample Pulse-Shifted P-M Generated Incident and Reflected 
Wave 




Figure 24. 



Sample Forward-Shifted P-M Generated Incident and 
Reflected Wave 




Figure 25. Original P-M Spectral Model 
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Figure 27. Log Amplitude of Pulse- Shif ted Combined P-M Generated Wave 
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Figure 28. Raw Phase of Pulse-Shifted Combined P-M Generated Wave 



AMPLITUDE * 

- 40.05 - 34.75 - 29.45 - 24.15 - 18.85 - 13.55 - 8.25 - 2.95 2.35 




FREQUENCY (CPS) 

Figure 29. Unwrapped Phase of Pulse-Shifted Combined P-M Generated Wave 
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Figure 30. Unwrapped Phase with Linear Trend Removed of Pulse-Shifted 
Combined P-M Generated Wave 
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Complex Cepstrum for Pulse-Shifted Combined P-M Generated 
Wave 



Figure 31. 




Figure 32. Power Cepstrum for Pulse-Shifted Combined P-M Generated Wave 
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Figure 34. Recovered Incident, True Incident and Combined Wave Time Series for Pulse-Shifted 
P-M Generated Wave 



AMPLITUDE - 

0.00 2.37 4.75 7.12 9.50 11.87 14.25 16.62 19.0 




FREQUENCY (CPS) 

Figure 35. Smoothed Spectrum of Pulse-Shifted Combined P-M Generated 
Wave 
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Figure 38. Smoothed Spectrum of Reflected P-M Generated Wave from 
Pulse-Shifted Combined P-M Generated Wave 
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Figure 39. Raw Spectrum of Forward-Shifted Combined P-M Generated Wave 
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Figure 40. Log Amplitude of Forward-Shifted Combined P-M Generated Wave 
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Figure 41. Raw Phase of Forward-Shifted Combined P-M Generated Wave 
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Unwrapped Phase of Forward-Shifted Combined P-M Generated 
Wave 



Figure 42. 
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Figure 43. Unwrapped Phase with Linear Trend Removed of Forward- 
Shifted Combined P-M Generated Wave 



AMPLITUDE - 

0.30 - 0.24 - 0.19 - 0.13 - 0.08 - 0.02 0.03 0.09 0.15 



86 




' o'. 00 

Figure 44. 
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Complex Cepstrum for Forward-Shif ted Combined P-M Generated 
Wave 
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Power Cepstrum for Forward-Shifted Combined P-M Generated 
Wave 



Figure 45. 
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Phase Cepstrum for Forward-Shifted Combined P-M Generated 
Wave 



Figure 46. 
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Figure 47* Recovered Incident, True Incident and Combined Wave Time Series for Forward- 
Shifted P-M Generated Wave 
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Figure 48. Smoothed Spectrum of Forward-Shifted Combined P-M Generated 
Wave 
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Figure 49. Smoothed Spectrum of Incident P-M Generated Wave 
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Figure 50. Smoothed Spectrum of Recovered Incident P-M Generated Wave 
from Forward-Shifted Combined P-M Generated Wave 
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Forward-Shifted Combined P-M Generated Wave 
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Figure 52. Smoothed Spectrum of Recovered Incident P-M Generated Wave 
from Forward-Shifted Combined P-M Generated Wave with 
Zeroes Added to End of Input Data 
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Figure 53. Smoothed Spectrum ot Reflected P-M Generated Wave from 
Forward-Shifted Combined P-M Generated Wave with Zeroes 
Added to End of Input Data 
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Fxgure 54. Complex Cepstrum for P-M Generated Wave, t Varying with 
Frequency, At=.333, a=.5 
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Figure 55. Filtered Complex Cepstrum for P-M Generated Wave, 
Varying with Frequency, At=.3333, a=.5 
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Figure 56. Complex Cepstrum for P-M Generated Wave, x Varying with 
Frequency, At=.7333, a=.5 
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Figure 57. Filtered Complex Cepstrum for P-M Generated Wave, x 
Varying with Frequency, At=.7333, ct=.5 
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Figure 58. Filtered Complex Cepstrum Using Filter Variation No. 1 on 
Complex Cepstrum Shown in Figure 54 
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Figure 59. Filtered Complex Cepstrum Using Filter Variation No. 2 on 
Complex Cepstrum Shown in Figure 54 
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Figure 60. Complex Cepstrum for P-M Generated Wave, t Varying with 
Frequency, a Varying from .4 to .6, At=.3333 
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Figure 61. Complex Cepstrum for P-M Generated Wave, x Varying with 
Frequency, a Varying from .7 to .9, At-. 3333 
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Figure 62. Sitka, Alaska, Breakwater Plan and Instrumentation Layout 
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Figure 63. Complex Cepstrum of 2048-Point Data Record of Combined 
Wave from Sitka Breakwater 
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Figure 64. Unwrapped Phase of 2048-Point Data Record of Combined Wave 
From Sitka Breakwater 
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Figure 65. Unwrapped Phase with Linear Trend Removed of 2048-Point 
Data Record of Combined Wave from Sitka Breakwater 
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Figure 66. Raw Spectrum of 1024-Point Data Record of Sitka Breakwater 
Combined Wave with 1024 Zeroes Added to End 
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Figure 67. 
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Log Amplitude of 1024-Point Data Record of Sitka 
Combined Wave with 1024 Zeroes Added to End 
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Figure 68. Unwrapped Phase of 1024-Point Data Record of Sitka Break- 
water Combined Wave with 1024 Zeroes Added to End 
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Figure 69. Unwrapped Phase with Linear Trend Removed of 1024-Point 
Data Record of Sitka Breakwater Combined Wave with 1024 
Zeroes Added to End 
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Figure 70. Complex Cepstrum of 1024-Point Data Record of Sitka 

Breakwater Combined Wave with 1024 Zeroes Added to End 
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Figure 71. Power Cepstrum of 1024-Point Data Record of Sitka Break- 
water Combined Wave with 1024 Zeroes Added to End 
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